Localization-delocalization transition of a polaron near an impurity 
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We solve the problem of polaron localization on an attractive impurity by means of direct-space 
Diagrammatic Monte Carlo implemented for the system in the thermodynamic limit. In particular 
we determine the ground state phase diagram in dependence on the electron-phonon coupling and 
impurity potential strength for the whole phonon frequency range. Including the quantum phonon 
dynamics we find and characterize a new phase which is missing in the zero phonon-frequency limit 
(adiabatic approximation), where self-trapped polarons are not localized at shallow impurities. We 
predict and show that in the vicinity of the localization transition a region with a mixture of weak- 
and strong-coupling spectral response is realized. 
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A general approach to the theoretical description of a 
particle in a bulk medium coupled both to bosonic ex- 
citations and the potential of imperfections is an impor- 
tant but notoriously hard problem that poses a real chal- 
lenge even to modern nonperturbative approaches [l[ . As 
yet only approximate results, relying, e.g., on dynamical- 
mean field theory exist. A central question in this con- 
text is the formation of three-dimensional (3D) polarons 
at impurities, or the Anderson localization of polarons in 
disordered media [1, 0, 0| ■ The overall importance of the 
physics of electron-phonon interaction in doped materials 
makes this issue of general interest for different areas of 
physics and technology. As a matter of fact the interplay 
between disorder and interaction effects is an important 
issue for contemporary materials design. For example 
high temperature superconductors [B|, Ig, 0] or materials 
with colossal magnetoresistance Q are doped Mott insu- 
lators where besides the coupling to bosonic excitations 
(phonons and magnons) disorder is present. 

In this Letter we present the exact solution to the po- 
laron problem in the presence of an attractive impurity 
in a 3D material. The accepted model for that situation 
is given by the Hamiltonian H = + ijW with 
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In H( ', U is the attractive impurity potential for the 
electron cj at site and b\ creates a dispersionless optical 
phonon with frequency w p h at Wannier i. describes 
the electron transfer oc t between nearest neighbor sites 
and local Holstein coupling to the phonons oc 7. 

In the absence of electron-phonon (el-ph) coupling 
(7 = 0), the critical U for particle localization at the 
impurity is t/ c (7 = 0) » 3.96 [§]; all energies are 



measured in units of t hereafter. In the adiabatic ap- 
proximation (AA), setting to p h = 0, the phase diagram 
in U — A coordinates, with the dimensionless coupling 
constant A = 7 2 /(6iw p h), was established in Ref. [loj ]. 
The phase boundary in AA, separating delocalized po- 
laron states from localized ones, crosses the [/-axis at 
U c {l = 0) and the A-axis at X C (U C — 0,w p h = 0) « 0.9. 
The latter crossing is a confusing property of the AA 
phase diagram since it implies that for el-ph couplings 
A > X C {U C = 0,Wph = 0) the polaron is localized even 
when U — 0. Quite the contrary, a particle is never 
localized in a translationally invariant lattice (U = 0) 
with quantum phonons (w p h > 0). Instead the particle 
undergoes only a crossover from the weak-coupling light 
polaron to a strong-coupling heavy polaron with small 
radius around a self-trapping coupling Ast U( • The AA 
erroneously equates Ast with the critical el-ph coupling 
strength A c required for polaron localization at U ^ 0. 
Therefore drastic differences between the exact result and 
that obtained in AA are expected, especially at small U c . 

Having this delicate situation in mind, we decided 
to study the full Hamiltonian {J)-© with quantum 
phonons. To this end we employ a new scheme combin- 
ing the Diagrammatic Monte Carlo (DMC) method in di- 
rect space [12J and the Stochastic Optimization Method 



for analytic continuation [12|, |13| which provides the 
approximation-free solution of the above problem with- 
out finite size errors and in zero temperature limit. Cal- 
culating the charge density distribution (CDD) around 
the impurity and the local density of states (LDOS) on 
the impurity site we establish the exact localization phase 
diagram for different phonon frequencies. We character- 
ize two novel polaronic regimes in a system with impuri- 
ties. The polaron at small U can be self-trapped though 
extended and not yet confined by shallow impurities. An- 
other regime, arising near the critical parameters for lo- 
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calization at the impurity, shows spectroscopic response 
like a mixture of spectra typical for weak, intermediate, 
and strong coupling. 

The direct space DMC method [14] can provide 
the direct space Green functions (GFs) in imaginary 
time (t) representation at zero temperature Gy(r) = 
(vac|cj (t)c; | vac) for the Hamiltonian (jT]I2J) by Feynman 
diagram expansion in the interaction representation 



Gij(r) = 



e -^ (0, TV 



(cj(r)cJ exp j- £H^( T ')dT' 



(3) 

The implementation of DMC [14| requires to keep in com- 
puter memory all GFs {Gy(r)} in direct space, which 
restricts the lattice to about 25 x 25 sites. To avoid this 
size limitation we calculate only quantities related to on- 
site GFs Gii(r). With our implementation of DMC we 
are able to calculate the on-site GFs at zero temperature 
for a 10 8 x 10 s x 10 s lattice, thereby avoiding any finite- 
size or finite temperature errors. A slight modification of 
Eq. Q, 



n(i) = 
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(ci(/3)4exp. 



H^(T')dT' 



introduces the estimator for the CDD at temperature 
T = /9 _1 . To make a calculation of the CDD feasible, we 
collect its statistics in a cube with 40 3 number of sites. 
Note that this strategy does not introduce finite-size er- 
rors because only the r = 0, (3 points of the partition 
function loop are confined to the 40 3 -cube while the dia- 
grams are free to sample all (10 8 ) 3 sites. 

The CDD estimator is effective for locating the local- 
ization parameters for large U only but, because of the 
requirement of finite temperatures, fails at small (7<i. 
Note that the path-integral quantum Monte-Carlo algo- 
rithm [l5j | , which is another method relevant for the prob- 
lem formulated above [3], has serious precision limits for 
the same reason. Hence, the only rigorous method to 
locate the localization point in the infinite system is to 
calculate the on-site zero temperature GF Gh(t), deter- 
mine the LDOS Li(w) = — 7r _1 ImGii(a;) by analytic con- 
tinuation 12, T3|, and check for the presence of a bound 
state in the LDOS L (oj) at the impurity site. 

To validate the new implementation of the DMC tech- 
nique, we located the critical U C (X = 0) ~ 3.96 by cal- 
culating the CDD, normalized to unity at the impurity 
site, around the impurity. It occurred that for U < U c 
the charge density does not decrease exponentially with 
distance from the impurity while for U > U c it does. Per- 
fect agreement is found between CDD obtained by DMC 
and that obtained in Ref. For U close to U c , however, 
determination of the LDOS Lq{lo) is a much more precise 
method, since the CDD requires finite temperatures. 

First let us demonstrate how trapped polaron states 
are determined using the LDOS. From the commutator 
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FIG. 1: (color online) Charge density at T = 0.001 (a) and 
LDOS at the impurity site for T — (c,d) for different values 
of U at ojpjj = 2 and A = 0.8. The arrow in panel (c) indicates 
the lower threshold (Th) of the spectrum at A = 0.8 and 
U — 0. Panel (b) shows LDOS at the impurity site at U = 2, 
cjph = 2 for A = (dashed line) and A = 0.8 (solid line). 
Statistical errorbars in panel (a) are less than 10~ 5 (i.e. much 
less than the point size). 



[H, c\] we find that, independent of the el-ph coupling A, 
the first moment Mi of the LDOS Li(uj) obeys 



/oo 
duj u>Li(uj) = -U6 it0 . 
-OO 



(5) 



In accordance with the sum rule ©, the LDOS at the 
impurity site shifts to lower energies with increasing U 
(Fig. fT};). The second moment Mi = f^°du> tu 2 Li(uj) in- 
creases with A (Fig. [TJa) and the overall LDOS broadens. 
Figure Q] shows the CDD (a) and LDOS (c,d) for various 
values of U at fixed A = 0.8 and uj p h — 2. Let us start the 
discussion with the case U = 0, where no localization is 
expected. We determine the lower border £>rh(A) of the 
LDOS Lq(uj) for given A (arrow in Fig. [TJc)). Increas- 
ing U the LDOS changes but there is no spectral density 
below ^Th(A) up to U < 1.90 (Fig. rjTc)). The CDD 
around the impurity, in accordance with the absence of 
a bound state in LDOS, does not show exponential de- 
crease too (Fig. Ha)). This gives another confirmation 
for the method employed here. In order to search for the 
localization-delocalization transition we proceed to larger 
values of U. For U > 1.95, the bound state appears 
below the threshold E Th (X) (Fig. [flc)), and the CDD 
decays exponentially (Fig. QJa)). In this way we obtain 
one transition point in the phase diagram (Fig. Ufa)), 
here U c {to ph = 2, A = 0.8) = 1.925 ± 0.025. Recall that, 
although the LDOS approach needs a very precise deter- 
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ruination (compare (Fig.QJc) and (d)), it is applicable for 
any values of A and U. On the contrary, the CDD method 
is fast but, due to requirement of finite temperature, not 
reliable at U <C 1. 
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FIG. 2: (color online) Phase boundaries between delocalized 
(left lower corner) and localized states of a polaron in A c — 
U c (a,b) and — U c (c) coordinates: adiabatic limit u>p^ — 
(thick solid line in panel (a)), Lj pll = 0.1 (squares), 0.5 
(triangles up), 1 (diamonds), 2 (triangles down), 4 (triangles 
left), 8 (crosses), 12 (stars), 16 (circles). The thick solid curve 
in panel (c) is obtained from Eq. ((6]). Solid, dashed and dotted 
lines in panel (b) are results obtained by the CBS method phi ] 
for a; p h = 8, 12, and 16. The values of A c are set exactly while 
the errorbars of the quantity U c are less than 3.0 x 10~ 2 . 



Next the phase diagram for polaron localization is pre- 
sented in Fig. [2] Using A c — U c coordinates (Fig. HJa), 
(b)), we see how our exact solution differs from the adi- 
abatic result (w p h — » 0, thick solid curve) for finite w p h- 
In Fig. (He), with g 2 — U c coordinates (g — j/u ph ), we 
show the deviation from the limiting phase boundary at 
cjph — > oo (thick solid curve), 
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FIG. 3: (color online) (a) Localization phase diagram and (b) 
average number of phonons JVav and the derivative dN av /dX 
for Lj ph = 4. (c) Dependence of the self-trapping coupling 
Ast on the phonon frequency UJph,- Errorbars are less than 
point size. 



We emphasize the excellent agreement of the novel 
DMC approach with the adiabatic limit [l(| and antia- 
diabatic limit Eq. © results, as well as with the data 
obtained by the Coherent Basis States (CBS) method 
[l6| for small U, which proves the validity of the imple- 
mentation. Note that the phase diagram obtained here 
is free from any substantial error and presents the first 
available solution of the Hamiltonian Eq. {T}-© for all 
parameter regimes. 
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U c (ujph = oo) = C/ c (7 = 0) exp(-gl). 



(6) 



(e) 













This relation is obtained by Lang-Firsov transformation, 
which renormalizes hopping, and accordingly the critical 
U c , as t — > £exp(— g 2 ). Note that a exponential relation 
between U c and g 2 (or A c = g 2 ui p ^/6t) is a characteristic 
property of the small [/ c < 1 and large A > 1 regime 
(Fig. H[b)), since for large cu p h, 

\c(U c , Lo ph ) ~ (io ph /6t) In [U e (n = 0)/U c ] + const . (7) 

It is indeed seen in Fig. Hfb) that the slope of the phase 
boundary increases with w p h- 



FIG. 4: (color online) LDOS for a localized polaron (a,b) at 
(7 = 2 and a delocalized polaron (c,d) at U = 1.9, where 
A = 0.8 and Wph = 2. Shown is the LDOS at the impurity 
site (0,0,0) (solid line), nearest neighbor (1,0,0) (dots), site 
(2,0,0) (dash-dots), (3,0,0) (dash-dot-dots), (4,0,0) (short- 
dashes), and for infinite distance from the impurity (short 
dots). The dashed line is the unperturbed LDOS (A = 0, 
U = 0). The inset in panel (b) gives the Z-factor of the LDOS 
5-peak of the bound state as a function of the distance to the 
impurity. Panel (e) shows the average number of phonons at 
(diamonds) and far from the impurity (squares) for u> p h = 1 
and U = 1.4. 
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Let us finally discuss two essential features of the phase 
diagram, which are entirely missing in the adiabatic ap- 
proximation. The first is realized at large el-ph couplings 
where the polarons are already self-trapped but not yet 
confined by shallow impurities (cross-hatched region in 
Fig. [21a)). The self-trapping coupling Ast, locating 
the crossover from the weak- to strong-coupling regime, 
can be defined, e.g., as the maximum of the derivative 
of the average number of phonons in the ground state 
iV av = ^&q = o^q=o^ with respect to the coupling con- 
stant A (Fig. El^b,c)). For small enough U <C f, and 
any given w p h, one finds a sector in the phase diagram 
(Fig- OH a)) where AsT(^ph) < A < A c (£7, cj p j,). This sec- 
tor defines the phase of self-trapped deconfined polarons, 
whose identification is obtained here for the first time. 

The second novel feature appears at moderate values 
of U close to the transition region between localized and 
extended states (line-shaded area in Fig.[3K). There, the 
spectral properties of a polaron are strongly position- 
dependent. In Fig. 0] we show the LDOS calculated at, 
and in the vicinity of the impurity. Both for a localized 
(a,b) and extended (c,d) polaron the LDOS at the impu- 
rity site strongly differs from that at the nearest neigh- 
bor site in the low-energy region (b,d). On the contrary, 
the overall features of LDOS at the nearest neighbor are 
very similar to those at infinite distance from the defect 
(Fig. [4^,c). This property points out how strongly the 
spectral properties depend on the value of the impurity 
potential at a given lattice site. Comparison of the av- 
erage number of phonons at the impurity with that in 
infinite distance to the impurity (Fig. |4j3) shows that for 
a wide range of parameters the lattice is weakly distorted 
far from the impurity while it is strongly deformed near 
the impurity. This demonstrates how impurities enhance 
the formation of small polarons. 

As a consequence it is expected that in a material 
with imperfections a mixture of behavior typical for 
weak-coupling polarons (far away from an impurity) and 
strong-coupling polarons (close to, or at, the impurity) 
occurs. Even though the impurity concentration can 
be small, the induced changes, e.g. in the spectral re- 
sponse, can be drastic. For example, for U = 1.90 and 
A = 0.8 (Fig. QJi) the charge density on the impurity site 
is four orders of magnitude larger than in the bulk of the 
system. Therefore, even a small impurity concentration 
n-i m 10 -3 suffices to entirely change the spectral prop- 
erties. For example, since photoemission [17J and optical 
conductivity [18j | spectra are a very different for weak and 
strong-coupling, one can expect very rich mixture of the 
spectral responses. 

In conclusion, introducing the direct space diagram- 
matic Monte Carlo in the thermodynamic limit we pre- 
sented the exact phase diagram for localization of a po- 
laron at an attractive impurity for all coupling strengths, 
values of the impurity potential, and phonon frequencies 
ranging from the adiabatic to the antiadiabatic regime. 



Most notably we characterize a novel phase where heavy 
polarons are mobile in the presence of shallow impuri- 
ties and predict complex spectral properties of the sys- 
tems close to the localization-delocalization transition. 
The present DMC method can be easily generalized to 
study more general situations, e.g. systems with long- 
range particle hopping, impurities with long-range at- 
tractive/repulsive potentials, or interfaces and layered 
structures, demonstrating the potential for future re- 
search. 
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